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Solubilization  of  the  armchair,  metallic  (10,10)  single-walled  carbon  nanotubes  (SWCNTs)  in 
toluene  is  modeled  using  molecular  dynamics  simulations.  Inter-  and  intra-molecular 
atomic  interactions  in  the  SWCNT  +  toluene  system  are  represented  using  COMPASS 
(Condensed-phased  Optimized  Molecular  Potential  for  Atomistic  Simulation  Studies),  the 
first  ab  initio  forcefield  that  enables  an  accurate  and  simultaneous  prediction  of  various 
gas-phase  and  condensed-phase  properties  of  organic  and  inorganic  materials. 

The  results  obtained  show  that  due  to  a  significant  drop  in  the  configurational  entropy  of 
toluene,  the  solvation  Gibbs  free  energy  for  these  nanotubes  in  toluene  is  small  but 
positive  suggesting  that  a  suspension  of  these  nanotubes  in  toluene  is  not  stable  and  that 
the  nanotubes  would  fall  out  of  the  solution.  This  prediction  is  consistent  with  experimental 
observations.  ©  2004  Kluwer  Academic  Publishers 


1.  Introduction 

Due  to  a  unique  combination  of  their  mechanical,  elec¬ 
trical  and  chemical  properties,  carbon  nanotubes  have 
been  investigated  very  aggressively  since  their  discov¬ 
ery  in  1991  [1].  Depending  on  the  fabrication  method 
used,  carbon  nanotubes  appear  either  predominantly  as 
single-walled  carbon  nanotubes  (SWCNTs)  or  as  multi- 
walled  carbon  nanotubes  (MWCNTs).  SWCNTs,  pre¬ 
dominantly  produced  in  carbon  ablation  and  arc  dis¬ 
charge  processes,  can  be  described  as  single  graphene 
sheets  rolled  up  into  a  cylinder  with  a  quasi-one- 
dimensional  crystal  structure.  Depending  on  their  diam¬ 
eter  and  the  spiral  conformation  (chirality),  SWCNTs 
can  be  either  semiconducting  or  metallic.  Mechani¬ 
cal  properties  of  SWCNTs  are  quite  remarkable;  their 
elastic  modulus  is  typically  above  lTPa,  and  they  can 
undergo  very  large  non-uniform  (even  highly  local¬ 
ized)  reversible  deformations.  Except  for  their  ends  and 
the  locations  of  topological  defects  (e.g.,  7-5-5-7  and 
Stone-Wales  defects),  SWCNTs  are  generally  not  very 
reactive. 

MWCNTs  are  generally  produced  during  thermal 
decomposition  of  carbon  precursors.  Due  to  a  weak 
inter-wall  bonding,  MWCNTs  have  generally  infe¬ 
rior  mechanical  properties  relative  to  those  of  the 
SWCNTs.  Their  electrical  properties  are  similar  to 
those  of  the  SWCNTs  although  they  can  not  be  eas¬ 
ily  correlated  with  their  chirality.  Chemical  properties 
of  the  MWCNTs  are  dominated  by  thestructure  of  their 


outer  wall  and  are,  hence,  similar  to  the  ones  of  the 
SWCNTs. 

While  carbon  nanotubes  have  been  perceived  as  hav¬ 
ing  a  great  potential  in  many  critical  applications  (e.g., 
field-emission  flat-panel  displays  [e.g.,  2],  novel  mi¬ 
croelectronic  devices  [e.g.,  3],  hydrogen  storage  de¬ 
vices  [e.g.,  4],  structural  reinforcement  agents  [e.g.,  5], 
and  chemical  and  electrochemical  sensors  [e.g.,  6]), 
the  lack  of  control  of  their  purity  and  diameter,  length 
and  chirality  distributions  during  fabrication  is  a  major 
current  obstacle  to  their  full  utilization.  During  fabri¬ 
cation  of  the  carbon  nanotubes,  other  materials  such  as 
carbon  onions  and  turbostratic/amorphous  graphite  are 
also  generally  produced.  Many  techniques  (e.g.,  oxi¬ 
dation)  have  been  proposed  for  separation  of  the  car¬ 
bon  nanotubes  from  the  unwanted  byproducts.  How¬ 
ever,  a  method  is  currently  lacking  for  separation 
of  the  carbon  nanotubes  according  to  their  diameter 
and/or  chirality,  the  key  geometrical  and  structural  pa¬ 
rameters  which  control  electronic  properties  of  these 
materials. 

Recently,  a  new  approach  for  carbon  nanotube  sepa¬ 
ration  and  purification  has  been  demonstrated  [7].  The 
method  is  based  on  non-covalent  sidewall  function¬ 
alization  (attachment  of  functional  groups/molecules 
on  the  outer  wall)  of  SWCNTs  and  MWCNTs 
with  it -conjugated  p o I y ( /; - p h c n I y c n e v i n y I c n e - co - 2 , 
5-dioctoxy-m-phenylenevinylene)  (PmPV-co-DOct- 
OPV)  polymer.  Pristine  carbon  nanotubes,  despite 
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their  hydrophobic  character,  are  found  not  to  be 
soluble  in  non-polar  solvents  such  as  toluene.  However 
after  functionalization  with  the  PmPV-co-DOctOPV, 
carbon  nanotubes  become  soluble  and  can  be  sus¬ 
pended  in  toluene.  While  solubilization  of  the  carbon 
nanotubes  through  sidewall  functionalization  can  be 
achieved  in  many  different  ways,  it  is  critical  that 
the  functionalization  is  of  a  non-covalent  character 
so  that  the  electronic  characteristics  of  the  carbon 
nanotubes  (governed  by  their  sp 2  hybridization)  are  not 
altered. 

In  a  series  of  papers,  Blau  and  co-workers  [7- 
16]  carried  out  a  detailed  experimental  and  computa¬ 
tional  investigation  of  SWCNT  and  MWCNT  solubi¬ 
lization  in  toluene  through  sidewall  functionalization 
with  PmPV-DOctOPV.  Their  atomistic  simulation  re¬ 
sults  showed  that  the  backbone  of  the  PmPV-DOctOPV 
molecule  adopts  a  relatively  flat  helical  structure  which 
is  governed  by  the  mefa-phenylene  linkage.  The  van 
der  Waals  interaction  between  octyloxy  groups,  on 
the  other  hand,  causes  these  groups  to  project  out¬ 
wards  from  the  helical  structure.  This  PmPV-DOctOPV 
conformation  is  shown  to  promote  adhesion  of  these 
molecules  to  the  nanotubes  and,  hence,  is  believed 
to  play  a  critical  role  in  the  nanotubes  solubilization 
process.  While  the  work  of  Blau  and  co-workers  [7- 
16]  has  resulted  in  a  major  improvement  of  our  un¬ 
derstanding  of  nanotubes  solubilization  with  sidewall 
functionalization,  two  important  points  were  not  ad¬ 
dressed  in  their  work:  (a)  why  hydrophobic  nanotubes 
are  not  soluble  in  nonpolar  solvents  like  toluene  and 
(b)  what  is  the  role  of  the  solvent  in  nanotubes  solu¬ 
bilization  by  sidewall  functionalization.  In  the  present 
work,  we  investigate,  using  atomistic  simulations,  the 
first  of  the  two  outstanding  issues,  i.e.,  solubilization 
of  the  carbon  nanotubes  in  toluene.  While  this  process 
can  be  expected  to  be  affected  by  the  nanotubes  chi¬ 
rality,  the  results  obtained  in  the  present  work  show 
that  the  chirality  effect  is  minor.  Hence,  only  the  re¬ 
sults  pertaining  to  solubilization  of  the  armchair  metal¬ 
lic,  (10,10)  SWCNTs  in  toluene  are  presented  in  this 
paper. 

The  organization  of  the  paper  is  as  follows:  A  de¬ 
scription  of  the  computational  cell,  the  computational 
method  and  the  inter-atomic  forcefield  potentials  used 
in  the  present  work  is  presented  in  Section  2.  The 
main  results  obtained  in  the  present  work  are  pre¬ 
sented  and  discussed  in  Section  3,  while  the  key  conclu¬ 
sions  resulted  from  the  present  study  are  summarized 
in  Section  4. 

2.  Computational  procedure 
2.1.  Computational  cell 

Computer  simulations  of  the  structure  of  pure  toluene 
are  carried  out  using  a  2.604  nm  by  2.604  nm  by 
2.604  nm  cubic  computational  cell  with  periodic 
boundary  conditions  applied  in  all  three  principal  direc¬ 
tions.  One  hundred  molecules  are  placed  in  the  cell  to 
obtain  the  average  density  of  toluene  essentially  equal 
to  its  experimental  counterpart  (867  kg/m3).  Toluene 
(C6H5CH3)  is  the  simplest  alkyl  benzene,  the  methyl- 
benzene,  in  which  one  of  the  hydrogen  atoms  in  the 


Figure  1  Geometrically  optimized  structure  of  a  toluene  molecule  and 
the  corresponding  atomic  charges. 

benzene  (CgHg)  is  replaced  by  a  methyl  group  (CH3). 
At  room  temperature,  toluene  is  in  liquid  state.  The 
geometrically  optimized  structure  of  toluene  along  with 
the  atomic  charges  is  shown  in  Fig.  1. 

To  model  the  interactions  between  isolated  (far 
spaced)  carbon  nanotubes  and  the  toluene  solvent,  a 
4.04  nm  by  4.04  nm  by  4.154  nm  rectangular  com¬ 
putational  cell  with  periodic  boundary  conditions  ap¬ 
plied  in  all  three  principal  directions  is  used.  A  single 
(10,10)  SWCNT  with  its  axis  aligned  in  the  "-direction 
is  placed  in  the  center  of  the  computational  cell.  Thus, 
the  SWCNT  is  considered  as  infinitely  long  and  the 
end-effects  associated  with  its  hemispherical  caps  are 
neglected.  The  end  effects  are  generally  considered  to 
play  a  minor  role  in  the  SWCNT  solubilization  process 
due  to  a  very  large  (ca.  100-1000)  length-to-diameter 
aspect  ratio  in  the  SWCNTs.  The  lengths  of  the  compu¬ 
tational  cell  in  the  x-  and  y-directions  are  sufficiently 
large  that  interactions  between  the  SWCNTs  in  the  ad¬ 
jacent  cells  can  be  neglected.  330  molecules  of  toluene 
are  placed  in  the  computational  cell  around  the  SWCNT 
and  the  size  of  the  computational  cell  in  the  x-  and 
y-directions  adjusted  until  the  toluene  density  in  the 
regions  far  away  from  the  SWCNT  reached  a  value 
close  to  the  average  experimental  density  of  toluene 
(867  kg/m3). 

2.2.  Computational  method 

The  (10,10)  SWCNT  +  toluene  system  discussed  in 
the  previous  section  is  modeled  using  classical  molecu¬ 
lar  dynamics  simulations  in  the  microcanonical  (NVE) 
ensemble.  To  ensure  stability  of  the  simulations  and 
energy  conservation,  a  constant  time  step  of  0.2  fs 
is  used  for  numerical  integration  of  the  equations  of 
motion. 

After  placing  the  nanotube  and  the  toluene  molecules 
into  the  computational  cell,  the  canonical  (NVT) 
molecular  simulations  are  first  carried  out  for  5  ps  un¬ 
til  the  system  is  equilibrated  at  a  desired  temperature 
of  300  K.  Temperature  control  is  achieved  by  velocity 
scaling  which  was  carried  out  every  1  fs.  Once  the  sys¬ 
tem  is  equilibrated,  microcanonical  (NVE)  simulations 
are  carried  out  for  additional  5  ps  and  the  data  colleted 
for  computation  of  the  system  properties. 
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Figure  2  A  snapshot  of  the  molecular  configurations  for:  (a)  pure  toluene  and  (b)  toluene  containing  SWCNTs. 


Snapshots  of  the  molecular  configurations  in  pure 
toluene  and  in  toluene  containing  a  SWCNT  is  shown 
in  Fig.  2a  and  b,  respectively. 

2.3.  Forcefield 

While  the  accurate  simulation  of  a  system  of  interacting 
particles  generally  entails  the  application  of  quantum 
mechanical  techniques,  such  techniques  are  computa¬ 
tionally  quite  expensive  and  are  usually  feasible  only 


in  systems  containing  up  to  few  hundreds  of  interacting 
particles.  In  addition,  the  main  goal  of  simulations  of 
the  systems  containing  a  large  number  of  particles  is 
generally  to  obtain  the  systems’  bulk  properties  which 
are  primarily  controlled  by  the  location  of  atomic  nu¬ 
clei  and  the  knowledge  of  the  electronic  structure,  pro¬ 
vided  by  the  quantum  mechanic  techniques,  is  not  crit¬ 
ical.  Under  these  circumstances,  a  good  insight  into 
the  behavior  of  a  system  can  be  obtained  if  a  reason¬ 
able,  physically-based  approximation  of  the  potential 


(forcefield)  in  which  atomic  nuclei  move  is  available. 
Such  a  forcefield  can  be  used  to  generate  a  set  of  system 
configurations  which  are  statistically  consistent  with  a 
fully  quantum  mechanical  description. 

As  stated  above,  a  crucial  point  in  the  atomistic  sim¬ 
ulations  of  multi-particle  systems  is  the  choice  of  the 
forcefields  which  describe,  in  an  approximate  manner, 
the  potential  energy  hypersuface  on  which  the  atomic 
nuclei  move.  In  other  words,  the  knowledge  of  force- 
fields  enables  determination  of  the  potential  energy  of 
a  system  in  a  given  configuration.  In  general,  the  poten¬ 
tial  energy  of  a  system  of  interacting  particles  can  be 
expressed  as  a  sum  of  the  valence  (or  bond),  faience* 
crossterm,  £Crossterm,  and  nonbond,  £n0nbond,  interaction 
energies  as: 


^angle-angle- torsion >  terms,  as. 

£crossterm  =  £bond-bond  T  £ angle-angle  T  £bond-angle 
T  foetid  bond-torsion  T  £middle  bond-torsion 
T  £angle- torsion  T  £angle-angle-tcrsion  (3) 

The  non-bond  interaction  term,  £nonbond:  accounts 
for  the  interactions  between  non-bonded  atoms  and  in¬ 
cludes  the  van  der  Waals  energy,  £vdWr  the  Coulomb 
electrostatic  energy,  £couiomb,  and  the  hydrogen  bond 
energy,  £hbond,  as: 

£non-bond  =  £vdW  T  ^Coulomb  T  £hbond  (4) 


£total  —  £valence  T  £crossterm  T  £non-bond  (1) 

The  valence  energy  generally  includes  a  bond  stretch¬ 
ing  term,  £bond>  a  two-bond  angle  term,  £angie,  a  dihe¬ 
dral  bond-torsion  term,  £torsion>  an  inversion  (or  an  out- 
of-plane  interaction)  term,  £oop,  and  a  Urey-Bradlay 
term  (involves  interactions  between  two  atoms  bonded 
to  a  common  atom),  £ub,  as: 


^valence 


=  £i 


bond 


+  E ; 


angle  ' 


+  £oop  +  £ub  (2) 


A  schematic  explaining  the  first  four  types  of  valence 
atomic  interactions  is  given  in  Fig.  3. 

The  crossterm  interacting  energy,  £Crossterm,  accounts 
for  the  effects  such  as  bond  lengths  and  angles  changes 
caused  by  the  surrounding  atoms  and  generally  in¬ 
cludes:  stretch-stretch  interactions  between  two  ad¬ 
jacent  bonds,  £bond-bond>  stretch-bend  interactions  be¬ 
tween  a  two-bond  angle  and  one  of  its  bonds,  £bond-angie> 
bend-bend  interactions  between  two  valence  angles 
associated  with  a  common  vertex  atom,  £angie-angie> 
stretch-torsion  interactions  between  a  dihedral  angle 
and  one  of  its  end  bonds,  Ecnll  bond-torsion>  stretch-torsion 
interactions  between  a  dihedral  angle  and  its  middle 
bond,  E middle _bond-torsion >  bend-torsion  interactions  be¬ 
tween  a  dihedral  angle  and  one  of  its  valence  an¬ 
gles,  E angle-torsion .  and  bend-bend-torsion  interactions 
between  a  dihedral  angle  and  its  two  valence  angles. 


L- 


b 
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Figure  3  A  schematic  of  the:  (a)  stretch,  (b)  angle,  (c)  torsion,  and 
(d)  inversion  valence  atomic  interactions. 


Inter-  and  intra-molecular  atomic  interactions  in  the 
SWCNT  +  toluene  system  described  in  the  previous 
section  are  modeled  using  COMPASS  (Condensed- 
phased  Optimized  Molecular  Potential  for  Atomistic 
Simulation  Studies),  the  first  ab  initio  forcefield  that 
enables  an  accurate  and  simultaneous  prediction  of  var¬ 
ious  gas-phase  and  condensed-phase  properties  of  or¬ 
ganic  and  inorganic  materials  [17-19].  The  COMPASS 
forcefield  uses  the  following  expression  for  various 
components  of  the  potential  energy: 


£bond  =  X)L K2(b  -  bof  +  K3(b  -  b0f 
b 


+  KA{b  -  b0)4i 

(5) 

E  angle  =  |tf2(0  ~  +  H3(9  -  00)3 

9 

+  H4(9  -  0O)4J 

(6) 

^torsion  =  ^2  [Vl  L1  “  C0S  (0  ~  0l)J 

0 

+  v2[l  -  cos  (20  -  0°) J 
+  V3|_1  -  cos  (30  —  0°)JJ 

(7) 

£ooP  =  K'/2 

X 

(8) 

£bond-bond  =  ^  ^  '  Efrfyib  b{))(b  Z?q) 

b  b' 

(9) 

£  angle-angle  =  ^  Fe0'{9  —  9o)(9'  —  9q ) 

9  9' 

(10) 

£bond-angle  =  ^  ^  '  £/;b (b  bo)(9  9q) 

b  9 

(11) 

^end_bond-torsion  —  ^  ^  ^o) [ ^ 1  COS  0 

b  0 

+  V2  cos  20  +  V3  cos  30] 

(12) 

^middle_bond-torsion  =  ^  ^  '  Fb'(j)(b  Z?q) 

b'  0 

x  [£1  cos  0  +  £2  cos  20  +  £3  cos  30] 

(13) 

£angle-torsion  =  ^  ^  '  £<90  (^  7p)[  V\  COS  0 

9  0 

+  V2  cos  20  +  V3  cos  30] 

(14) 
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angle-angle-torsion  — 


J2 12  K* eei 

tpee' 


x  cos  (pie-  9o)(Q'-0'0)  (15) 


\  '  <M]_ 

^Coulomb  —  /  _ 

Z—  er\ i 


1  >j 


-vdW 


-E 


5 


<>;  Lrij 


_u 

,-6 

'>.i  J 


(16) 


(17) 


where  Z?  and  //  are  the  bond  lengths,  0  the  two-bond 
angle,  0  the  dihedral  torsion  angle,  /  the  out  of  plane 
angle,  q  the  atomic  charge,  s  the  dielectric  constant,  ry 
the  i- j  atomic  separation  distance,  bo,  Kt(i  =  2-4),  0q, 
Hi  (i  =  2-4),  0?  (i  =  1-3),  Vi (i  =  1-3),  FbW,  b'0.  Foe s 
$0’  ^b'0 ^  F{(i  —  1—3,  Fq^,  K^qq',  Ay,  and 

are  the  system  dependent  parameters  implemented  into 
Discover  [20],  the  atomic  simulation  program  used  in 
the  present  work. 


3.  Results  and  discussion 
3.1.  Pure  liquid  toluene 

To  analyze  the  structure  of  the  pure  toluene  (toluene  in 
the  absence  of  a  SWCNT),  a  radial  distribution  func¬ 
tion  and  several  angular  parameters  plots  are  gener¬ 
ated  using  the  simulation  results  and  analyzed  in  this 
section. 

A  number  density  radial  distribution  function  for 
pure  toluene  is  displayed  in  Fig.  4a.  The  radial  posi¬ 
tions  indicated  are  associated  with  the  center  of  mass 
of  the  toluene  molecules.  A  simple  observation  of  the 
results  shown  in  Fig.  4a  indicates  the  presence  of  a 
well-defined  first  solvation  shell  and  a  diffuse  second 
solvation  shell  surrounding  toluene  molecules.  The  re¬ 
sults  obtained  are  very  similar  to  the  ones  reported  by 
Kim  and  Lee  [21]. 

A  probability  density  plot  for  the  angle  the  bond  be¬ 
tween  the  methyl  group  and  the  benzene  ring  makes 
with  the  benzene-ring  plane  is  shown  in  Fig.  5a.  It  is 
seen  that  while  the  most  probable  value  of  this  angle 
is  zero  deg.  (the  value  in  the  geometrically  optimized 
structure  of  an  isolated  toluene  molecule),  a  substantial 
fraction  of  toluene  molecules  have  this  angle  as  large 
as  10  deg.  This  finding  suggests  that  interactions  be¬ 
tween  the  toluene  molecules  cause  minor  distortions  in 
the  structure  of  toluene  molecules. 

A  probability  density  plot  for  the  angle  the  benzene¬ 
ring  normal  makes  with  one  of  the  principal  coordinates 
(arbitrarily  chosen  as  the  --axis)  is  plotted  in  Fig.  5b. 
A  dashed  line  is  also  shown  in  Fig.  5b  which  rep¬ 
resents  the  corresponding  radial  distribution  function 
(  f(6)  =  7r/360sin($),  where  6  is  in  degrees)  associ¬ 
ated  with  a  random  distribution  of  toluene  molecules. 
By  comparing  the  simulation  results  with  the  dashed 
curve,  it  can  be  concluded  that  the  structure  of  pure 
toluene,  past  the  first  solvation  shell,  is  quite  random. 

A  probability  density  plot  for  the  angle  the  bond  be¬ 
tween  the  methyl  group  and  the  benzene  ring  makes 
with  the  z-axis  is  plotted  in  Fig.  5c.  Again,  the  cor¬ 
responding  probability  density  function  for  a  random 
distribution  of  toluene  molecules  is  displayed.  It  is  seen 


(a) 


(b) 

Figure  4  Number-density  radial  distribution  function  for:  (a)  pure 
toluene  and  (b)  toluene  containing  a  SWCNT. 


that  with  respect  to  this  angular  parameter,  the  structure 
of  the  pure  toluene  appears  also  quite  random. 

3.2.  Single  SWCNT  solubilized  in  toluene 
3.2. 7.  The  structure  of  toluene  in  the 
presence  of  a  SWCNT 

The  structure  of  toluene  in  the  presence  of  a  SWCNT 
is  analyzed  using  the  number-density  radial  distribu¬ 
tion  function  for  carbon  and  hydrogen  atoms  in  the 
toluene,  as  well  as  by  examining  the  orientation  of 
the  toluene  molecules  with  respect  to  the  axis  of  the 
nanotube.  The  effect  of  SWCNT  on  the  structure  of 
individual  toluene  molecules  (primarily  the  orientation 
of  the  methyl/benzene-ring  bond)  is  also  studied. 

The  number-density  radial  distribution  function  re¬ 
sults  for  toluene  in  the  presence  of  a  SWCNT  is  shown 
in  Fig.  4b.  Separate  plots  are  generated  for  carbon  and 
hydrogen  in  the  benzene  ring  and  the  carbon  and  hy¬ 
drogen  in  the  methyl  group.  The  results  displayed  in 
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Figure  5  Angular  parameters  characterizing  the  structure  of  pure  toluene.  Please  see  text  for  details. 


Fig.  4b  can  be  summarized  as  follows: 

(a)  Toluene  molecules  surrounding  the  SWCNT 
form  a  well-defined  first  and  a  more-diffuse  second  sol¬ 
vation  shell;  and 

(b)  The  orientation  of  toluene  molecules  in  the  first 
solvation  shell  is  such  that  the  benzene-ring  portion  of 
these  molecules  is  closer  to  the  nanotube  axis  then  the 
methyl-group  part  of  these  molecules.  A  more  complete 
picture  about  the  orientation  of  toluene  molecules  is 
obtained  using  the  analysis  presented  below. 

A  probability  density  plot  for  the  angle,  the  methyl 
bond  makes  with  the  benzene  plane  for  toluene 
molecules  in  the  first  solvation  shell  is  displayed  in 
Fig.  6a.  By  comparing  the  results  displayed  in  Figs  6a 
and  5a,  it  can  be  concluded  that  the  presence  of  SWCNT 
does  not  significantly  affect  the  toluene  intra-molecular 
angle  in  question,  and  hence,  most  likely  does  not  sig¬ 
nificantly  affect  the  structure  of  toluene  molecules  rel¬ 
ative  to  that  in  the  pure  solvent. 

A  probability  density  plot  for  the  angle  between  the 
benzene  plane  normal  of  the  toluene  molecules  in  the 


first  solvation  shell  and  the  carbon  nanotube  axis  is 
displayed  in  Fig.  6b.  The  results  displayed  suggest  that 
the  toluene  molecules  in  the  first  solvation  shell  tend  to 
align  themselves  in  such  a  way  that  their  benzene  plane 
is  parallel  to  the  SWCNT  axis.  When  an  analogous  plot 
was  generated  for  the  second  solvation  shell  (the  results 
not  shown  for  the  brevity),  no  significant  reorientation 
of  the  toluene  molecules  relative  to  the  structure  of  pure 
toluene  is  observed. 

A  probability  density  plot  for  the  angle,  the  bond 
between  the  methyl  group  and  the  benzene  ring  makes 
with  the  SWCNT  axis,  for  the  toluene  molecules  in 
the  first  solvation  shell  is  displayed  in  Fig.  6c.  It  is 
seen  that  the  most  probable  values  for  this  angle  are 
~30  deg.  and  120  deg.  This  finding  suggests  that  one  of 
the  carbon-carbon  bonds  in  the  benzene  ring  involving 
the  carbon  atom  of  the  benzene  ring  which  is  bonded  to 
the  methyl  group,  is  nearly  orthogonal  to  the  SWCNT 
axis. 

A  probability  density  plot  for  the  angle  between  ra¬ 
dial  vectors  connecting  the  SWCNT  axis  and  the  cen¬ 
ter  of  the  benzene  ring  and  the  benzene  ring  normal, 
for  the  toluene  molecules  in  the  first  solvation  shell,  is 
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Figure  6  Angular  parameter  characterizing  the  structure  of  toluene  in 


Angle  Between  Benzine  Plane  Normal  and  Radial  Direction 

presence  of  a  SWCNT.  Please  see  text  for  details. 


displayed  in  Fig.  6d.  The  results  displayed  in  this  figure 
combined  with  the  ones  shown  in  Fig.  6b  suggest  that 
the  benzene  plane  of  the  toluene  molecules  in  the  first 
solvation  shell  is  tangential  to  the  cylindrical  walls  of 
the  first  solvation  shell. 


3.2.2.  SWCNT  solvation  energy  in  toluene 

Molecular  dynamic  simulation  results  averaged  over 
the  last  5  ps  of  the  simulation  time  are  also  used  to  an¬ 
alyze  the  energetics  of  the  process  of  the  introduction 
of  a  SWCNT  from  the  vacuum  into  the  bulk  toluene 
solvent.  The  total  energy  change  accompanying  this 
process,  A£t+swcnt,  generally  referred  to  as  the  sol¬ 
vation  energy,  has  three  contributions,  i.e. : 

A£t+swcnt  =  A£t  +  A£Swcnt  +  ^t+cnt  (18) 

where,  A  Ej  represents  a  change  in  the  solvent  (toluene) 
energy  accompanying  the  creation  of  a  cylindrical  cav¬ 
ity  needed  to  accommodate  the  SWCNT.  A  Ej  (ana¬ 


lyzed  below  in  more  details)  has  two  contributions  it¬ 
self,  one  associated  with  a  change  in  the  bulk  struc¬ 
ture  of  the  solvent  caused  by  the  presence  of  the  solute 
(SWCNT)  and  the  other  associated  with  the  surface  en¬ 
ergy  of  the  cylindrical  cavity.  A  f  swcnt  in  Equation  18 
represents  a  change  in  the  SWCNT  energy  due  to  sol¬ 
vent/solute  interactions.  A  £ swcnt  is  primarily  the  re¬ 
sult  of  the  changes  in  the  carbon-carbon  inter-atomic 
distances  within  the  SWCNT  caused  by  interactions 
between  the  SWCNT  and  the  toluene  molecules  in  the 
first  solvation  shell.  £ swcnt  is  the  energy  of  the  inter¬ 
action  between  the  atoms/molecules  of  the  solvent  and 
the  solute. 

Based  on  the  atomic-scale  simulations  results  ob¬ 
tained,  the  following  energy  change  accompanying 
introduction  of  a  SWCNT  in  toluene  is  obtained: 

A£t+swcnt  =  Et+swcnt  —  Ej  —  E  swcnt 
=  128,805.0  -  (-14,472.8) 

-  145,592.7  =  -2,314.9  kJ/mol 
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where  Ej  and  £SWcnt  are  the  molar  energies  of  the 
bulk  toluene  and  the  pristine  SWCNT  in  the  vacuum, 
respectively.  Superscript  zero  is  used  in  the  remainder 
of  the  paper  to  denote  various  quantities  of  the  pure 
substances  in  vacuum.  It  should  be  noted  that  the  term 
mole  used  in  the  present  paper  pertains  to  the  Avo- 
gadro’s  number  of  computational  cells.  In  order  words, 
all  species  within  a  computational  cell  are  treated  as 
one  (complex)  molecule. 

The  average  changes  in  the  molar  energies  of  the  bulk 
toluene  and  the  SWCNT  are  found  as  follows: 

A Et  =  Ej  -  Ej  =  -14,078.2  -  (-14,472.8) 
=  394.6  kJ/mol,  and 

A£sWCNT  =  £swcnt  —  ^SWCNT 

=  145,596.9  -  145,992.7  =  4.2  kJ/mol. 

The  interaction  energy  between  the  toluene  and  the 
SWCNT,  £swcnt’  's  computed  as: 


where  the  energies  are  expressed  in  J/mol.  NA  is  the 
Avogadro’s  number  and  A  the  exposed  surface  area  in 
m2.  The  value  of  the  surface  energy  (0.016  J/m2)  is  not 
in  very  good  agreement  with  its  experiment  counterpart, 
0.018  J/m2  [22], 

Since  the  cavity  that  accommodates  a  SWCNT 
within  the  solvent  is  cylindrical,  the  curvature  of  its 
surface  can,  in  principle,  have  an  effect  on  the  surface 
energy.  To  quantify  this  effect,  attempts  were  made  to 
carry  out  geometrical  optimization  of  the  toluene  struc¬ 
ture  containing  a  cylindrical  cavity.  Unfortunately,  the 
attempts  were  not  successful  since  the  cavity  was  un¬ 
stable  and  was  quickly  filled  by  the  surrounding  toluene 
molecules.  Consequently,  a  different  strategy  was  pur¬ 
sued.  The  effect  of  the  surface  curvature  on  the  surface 
energy  was  assessed  by  comparing  the  surface  energies 
of  unrelaxed  (geometrically  un-optimized)  structures 
of  the  slab  toluene  and  the  one  containing  a  cylindrical 
cavity.  The  average  (unrelaxed)  surface  energy  of  a  flat 
toluene  surface  is  calculated  as: 


Ct+swcnt  =  A^t+swcnt  —  A£t  —  A£swcnt 
=  -2,314.9-394.6-4.2 
=  -2,713.7  kJ/mol 

Since  the  total  energy  change  associated  with 
the  introduction  of  a  SWCNT  into  the  toluene 
(A  Cj+swcnt  =  —2,314.9  kJ/mol)  is  negative,  our 
preliminary  finding  is  that  the  SWCNTs  should  be  solu¬ 
ble  in  this  solvent.  However,  to  obtain  a  definite  answer 
to  this  question,  one  must  also  consider  the  entropy 
change  accompanying  the  introduction  of  a  SWCNT 
into  the  toluene.  An  analysis  of  the  SWCNT  solvation 
entropy  in  toluene  is  next  section. 

The  results  presented  above  show  that  the  change  in 
the  toluene  energy  due  to  introduction  of  a  SWCNT, 
A  Ej  =  394.6  kJ/mol,  makes  a  significant  contribu¬ 
tions  to  the  SWCNT  solvation  energy.  As  mentioned 
above,  this  component  of  the  solvation  energy  arises 
from  the  creation  of  a  surface  in  the  cylindrical  cav¬ 
ity  and  due  to  the  rearrangement  of  toluene  molecules 
(primarily  the  ones  in  the  first  and  the  second  solvation 
shells  surrounding  the  SWCNT).  To  assess  the  relative 
magnitudes  of  the  two  contributions  to  A  Ej.  the  sur¬ 
face  energy  of  toluene  is  first  computed  by  carrying  out 
molecular  dynamics  simulations  described  earlier  but 
for  a  slab  (instead  for  a  bulk)  geometry  of  the  toluene 
computational  cell.  The  slab  is  obtained  by  removing 
the  periodic  boundary  conditions  in  one  of  the  prin¬ 
cipal  directions  of  the  bulk  computational  cell  which 
results  in  the  formation  of  two  parallel  surfaces  per 
computational  cell.  The  surface  energy  of  toluene  is 
next  obtained  as: 


£  surf,  flat 


(E 


o,slab 

T 


Ej) 


NaA 


-4335.1  x  103  +  4468.4  x  103 
NA  x  2  x  (26.04  x  10-10)2 

0.016  J/m2 


E 


surf,  flat,  unrel 
T 


(£^'slab,  unrel  -  Ej) 

N^A 

757.4  x  103  -  88.8  x  103 
Na  x  2  x  (26.04  x  10-10)2 

0.082  J/  m2 


while  the  corresponding  surface  energy  of  the  cylindri¬ 
cal  cavity  is  computed  as: 


E 


surf,  cyl,  unrel 
T 


(£°’cyl  -  Ej) 

NaA 

4330.9  x  103  -  3612.8  x  103 
NAx2n  x  9.9  x  26.04  x  lO”20 

0.074  J/m2 


The  two  surface  energy  values,  which  differ  by  about 
around  10%,  suggest  that  the  curvature  most  likely  does 
not  affect  the  surface  energy  (in  the  relaxed  structures) 
very  significantly. 

The  energy  change  of  the  toluene  arising  from  for¬ 
mation  of  a  cylindrical  cavity  which  accommodates  the 
SWCNT  is  next  calculate  by  multiplying  the  calculated 
surface  energy  with  the  surface  area  of  the  cylindrical 
cavity  as: 


£™rf'cyl  =  0.016  x  Na  x  (2tt  x  9.9  x  26.049 
x  10“20)  x  10“3  =  254.1  kJ/mol 


A  comparison  of  this  value  with  A  Ej  =  394.6  kJ/mol 
suggests  that  the  energy  change  of  the  toluene  associ¬ 
ated  with  formation  of  a  cylindrical  cavity  represents  al¬ 
most  two  thirds  of  the  total  energy  change  of  the  toluene 
asoociated  with  the  introduction  of  a  SWCNT. 
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3.2.3.  SWCNT  solvation  entropy  in  toluene 

The  molar  entropy  of  pure  toluene  can  be  expressed  as  a 
sum  of  the  corresponding  ideal-gas  molar  entropy,  Sj, 
and  the  molar  excess  entropy,  S™,  as: 

ST  =  48  +  5tx  (19) 

where  all  three  quantities  in  Equation  19  are  evaluated 
at  a  given  volume,  V,  and  temperature,  T,  conditions 
of  the  system. 

The  ideal-gas  molar  entropy  can  be  expressed  as: 

4  =  f-,uln(^)  (20) 

where  Ru  is  the  universal  gas  constant,  f>j  the  number 
density  of  toluene  and  Ax  the  corresponding  thermal 
de  Broglie  wavelength  defined  as: 

<2i» 

where  h  is  the  Planck’s  constant,  mj  is  the  molecular 
mass  of  toluene,  and  A'b  the  Boltzmann’s  constant. 

The  molar  excess  entropy  can  be  defined  as  a  sum 
of  two-,  three-,  and  higher  order  many-body  contribu¬ 
tions.  Computational  studies  of  liquid  metals  [23],  no¬ 
ble  gases  [24] ,  hard-sphere  fluids  [24, 25  ] ,  and  Lennard- 
Jones  fluids  [25,  26]  at  liquid-like  densities,  showed 
that  the  molar  excess  entropy  is  dominated,  and  hence 
can  be  reasonably  well  approximated,  by  the  two-body 
term.  This  simplification  is  utilized  in  the  present  work. 

If  the  toluene  molecules  are  considered  tentatively 
as  spherically  shaped  for  simplicity,  the  toluene  two- 
body  (pair)  correlation  function  depends  only  on  the 
relative  position  vector,  r,  of  the  molecules.  Under  such 
approximation,  the  molar  excess  entropy  of  toluene  can 
be  expressed  as: 

5ex  =  _^T  j  g(r)ln(£(r))cb 

+  Udv  (22) 


where  g(r )  is  the  pair  correlation  function  and  the  inte¬ 
gration  is  carried  out  over  the  system  volume. 

Under  a  similar  assumption,  the  molar  entropy  of  a 
solution  consisting  of  SWCNTs  suspended  in  toluene 
can  be  expressed  as: 

St+swcnt 

=  - 2 - ^u(d  —  -^swcNTllniPTAj) 

5i?uVswCNT  D  ,  /  .  3  \ 

H - 2 - «u-VSWCNTln[PsWCNTASWCNTj 

*u(l  —  *swcnt)Pt 
2 

X  /gx,(r)ln[gT,T(r)]dn-/ [gT,T(r)  -  l]dnj 


/ 


—  /^u-VSWCNTPT 

[gswcNT,T(r )  —  1]  dn 


/' 

1/ 


gswcNT.TU )  lnlgswcNT.TU )]  du 
^u-rsWCNTPSWCNT 


x  |  /  gswcNT,swcNT(r )  lnfgswcNT.swcNT^)]  dv 


x  J  [£swcnt,swcntU )  —  1]  dr 


(23) 


where  a'swcnt  is  the  mole  fraction  of  SWCNTs  in 
the  solution,  and  gx,T,  gT.  swcnt,  and  gswcNT.swcNT 
are  the  toluene-toluene,  toluene-SWCNT  and  SWCNT- 
SWCNT  pair  correlation  functions,  respectively. 

The  partial  molar  entropy  of  the  solute,  5'swcnt-  is 
obtained  by  differentiating  Equation  23  with  respect 
to  IVswcnt  =  MwV'swcnt  under  constant  T,  V  and 
Ay  =  Aav(1  —  vswcnt)  where  Aav  is  the  Avogadro’s 
number.  In  the  limit  of  infinite  dilution,  this  procedure 
yields: 


SsWCNT 

=  -jp-  ~  -^uln(psWCNT^'-sWCNT. 


^  swcnt) 

J  gSWCNT.Tln  gSWCNT.Tdn 

^  LI  Pf 


—  ^uPT 
x  —  J  (gSWCNT.T  —  1)  du 


f  /  dgT.T  \ 
J  V  fiPSWCNT  / 


ln(gx,T)  dn  (24) 


t.v,nt 


The  first  two  terms  in  Equation  24,  reflect  the  ideal- 
gas  contributions,  while  the  last  two  terms  are  associ¬ 
ated  with  the  solute- solvent  correlation  function.  Fol¬ 
lowing  Ashbaugh  and  Paulaitis  [27],  the  last  term  in 
Equation  24  can  be  neglected  since  it  generally  makes 
a  small  contribution  to  the  solvation  entropy. 

Also,  following  the  procedure  described  by  Ash¬ 
baugh  and  Paulaitis  [27],  the  standard  solvation  entropy 
ASt+swcnt  is  related  to  the  partial  molar  entropy  of 
solute,  SswcnTl  as: 


ASt+SWCNT  =  SsWCNT  —  Ru  ln(psWCNTAsWCNT) 


-  RuTotj 


3  Ru 
2 


(25) 


where  ax  =  9.6  x  10-(7K  is  the  thermal  expansion 
coefficient  of  toluene.  Substitution  of  Equation  24  into 
Equation  25  yields: 


ASx+swcnt 


=  —RuTo/j 


RuPt 


J  gSWCNT.Thl  gSWCNT.T  dl> 


J  (gSWCNT.T  —  1)  du 


(26) 


The  second  term  in  the  Equation  26,  is  generally  re¬ 
ferred  to  as  the  solute-solvent  correlation  entropy  and 
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under  the  assumption  of  spherical  shapes  of  the  sol¬ 
vent  and  solute  molecules,  the  pair  correlation  function 
depends  only  on  the  radial  distance,  r.  However,  in  the 
present  case  where  toluene  molecules  are  plate-like  and 
the  SWCNT  has  a  cylindrical  shape,  the  pair-correlation 
function  involves  additional  degrees  of  freedom  which 
are  required  to  describe  solute-solvent  orientations  and 
conformations.  In  the  present  work,  it  is  assumed  that 
within  the  first  solvation  shell  of  the  SWCNT,  the 
solute-solvent  pair  correlation  function  can  be  factor¬ 
ized  into  a  radial  part  and  an  orientation-dependent  part. 
Beyond  the  first  solvation  shell,  on  the  other  hand,  the 
pair  correlation  function  is  assumed  to  involve  only  the 
radial  part.  Thus,  the  toluene-SWCNT  pair  correlation 
function  is  defined  as: 

£t,swcnt(t  ct,  p,  y) 

__  |^T,SWCNT(r)^T,lwCNT(“’  P’  T )  r  -  rsh 
I^tIsWCNtW  r  >  rsh 

where  a,  ft  and  y  are  Euler  angles  (defined  below),  rsh 
the  radius  of  the  first  solvation  shell,  while  superscripts 
rad  and  ang  are  used  to  denote  the  radial  and  the  an¬ 
gular  (orientation-dependent)  components  of  the  pair 
correlation  function.  For  the  toluene-SWCNT  pair  cor¬ 
relation  function  given  by  Equation  27,  the  solvation 
entropy  originally  defined  by  Equation  26  becomes: 

ASt+swcnt 

=  —  KPl  |  j  k+SWCNT^' ) 

X  ^[StJSWCNtW]  df  ~  J  [<?Ta<SWCNT(r)  —  l]  du  J 

-  Vsh  jjj  StIwCNT^  P-  Y) 

x  ln[gTnswcNT(a’  P>  y)]  sin/d  da  x  dp  x  dyj 

(28) 

where  /\x<,nn  is  the  normalization  factor  for  the  angular 
part  of  the  pair-correlation  function  (defined  below)  and 
the  Vsh  term  represents  the  integral  of  the  radial  part  of 
the  pair  distribution  function  over  the  volume  of  the 
first  solvation  shell. 

The  orientation  of  each  toluene  molecule  with  re¬ 
spect  to  the  nanotube  axis  is  defined  using  three  Eu¬ 
ler  angles.  Before  quantifying  Euler  angles  for  each 
toluene  molecule,  a  Cartesian  coordinate  system  is  at¬ 
tached  to  such  a  molecule  as  shown  in  Fig.  7.  The  x- 
and  z-axes  of  the  coordinate  system  are  aligned  with 
the  methyl-bond  direction,  Rm,  and  the  benzene-plane 
normal,  Rn,  respectively.  The  y-axis  in  the  same  coor¬ 
dinate  system  is  then  obtained  as  y  =  Ra  x  Rm,  where 
x  is  used  to  denote  a  cross  product  of  two  vectors.  It 
should  be  noted  that  since  toluene  molecules  are  dis¬ 
torted,  as  shown  in  Section  3.2.1,  a  simple  procedure 
was  developed  to  determine  the  “effective”  benzene- 
plane  normal  and  the  “effective”  methyl-bond  direction 


Figure  7  The  SWCNT-based  X-Y-Z  and  the  toluene  molecule-based 
x-y-z  Cartesian  coordinate  systems  used  in  the  present  work. 

for  toluene  molecules  before  the  coordinate  system  de¬ 
scribed  above  can  be  constructed.  Next,  a  Cartesian 
coordinate  system  associated  with  the  SWCNT  is  next 
constructed  in  the  following  way:  the  X-axis  is  aligned 
with  the  SWCNT  axis,  while  the  Z-axis  is  aligned  with 
a  radial  direction  originating  from  the  SWCNT  axis  and 
passing  through  the  center  of  the  benzene  plane  of  the 
toluene  molecule  in  question.  The  y-axis  is  then  de¬ 
fined  as  Y  =  Z  x  X.  Once  the  two  coordinate  systems 
are  established,  the  three  Euler  angles  are  defined  as 
following:  a  is  defined  as  the  angle  between  zxZ  and 
y,  p  as  the  angle  between  z  and  Z  and  y  as  the  angle 
between  z  x  Z  and  Y.  The  ranges  of  the  three  Euler 
angles  are  0  <  a  <  2n,  0  <  p  <  n  and  0  <  y  <  2n. 
Since  a  differential  element  of  the  solid  angle  is  de¬ 
fined  as:  d£2  =  sin  p  ■  da  ■  d/d  •  dy,  and  the  integration 
of  dfi  between  the  limits  of  a,  p  and  y  stated  above 
yields  A xoml  =  8tt2,  A xomi  is  the  appropriate  normal¬ 
ization  factor  in  the  angular  part  of  the  pair-correlation 
function. 

The  solvation  entropy  for  SWCNTs  in  toluene  is  then 
calculated  using  Equation  28  and  the  atomistic  simu¬ 
lation  results  for  the  last  5  ps  of  the  simulation  time. 
Integrals  appearing  in  Equation  28  are  evaluated  nu¬ 
merically  using  the  trapezoid  rule.  Due  to  the  cylin¬ 
drical  geometry  of  the  SWCNT,  the  volume  element 
appear  in  Equation  28  is  defined  as  dv  =  2nazrdr, 
where  az  is  the  computational-cell  lattice  parameter 
in  ^-direction.  The  radial  and  the  angular  components 
of  the  toluene-SWCNT  pair-correlation  function  are 
evaluated  numerically  by  constructing  the  correspond¬ 
ing  frequency  functions.  The  results  of  this  procedure 
yielded  the  solvation  entropy  for  SWCNT  in  toluene 
of  ASt+swcnt  =  —8.6  kJ/mol/K.  Using  this  value  of 
solvation  entropy,  the  solvation  Gibbs  free  energy  is 
evaluated  as: 

AGt+swcnt  =  A£t+Swcnt  —  T  ASt+swcnt 
=  -2042.0  -  298.0  x  (-8.6) 

=  495  kJ/mol 

where  one  mole  is  defined  as  the  Avogadro’s  number 
of  computational  cells.  Since  this  value  of  the  solva¬ 
tion  Gibbs  free  energy  is  small  (one  computational  cell 
contains  330  toluene  molecules)  but  positive,  one  would 
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expect  that  SWCNTs  could  be  relatively  easily  suspend 
in  toluene  but  the  resulting  suspension  would  not  be 
stable  and,  given  enough  time,  SWCNTs  would  drop 
out  of  solution.  This  prediction  is  fully  consistent  with 
the  experimental  observation  made  by  in  het  Panhuis 
et  al.  [28]  during  their  attempts  to  solubilize  SWCNTs 
in  toluene. 

4.  Conclusions 

Based  on  the  results  obtained  in  the  present  work,  the 
following  main  conclusions  can  be  drawn: 

1.  Introduction  of  the  SWCNTs  into  toluene  is  as¬ 
sociated  with  a  relatively  small  but  negative  solvation 
energy.  The  largest  contribution  to  be  solvation  energy 
arises  from  the  SWCNT-toluene  interactions. 

2.  Interactions  between  SWCNTs  and  toluene 
molecules  result  in  major  reorganization  of  the  toluene 
molecules  and  the  associated  conformation  gives  rise 
to  a  substantial  decrease  in  the  configurational  entropy 
of  the  toluene. 

3.  The  solvation  Gibbs  free  energy  of  SWCNTs  in 
toluene  at  room  temperature  is  relatively  small  but 
positive  suggesting  that  suspension  of  SWCNTs  in 
toluene  is  not  stable.  This  prediction  is  consistent  with 
experimental  observations. 
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